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ABSTRACT 

Aims. We probe the effect of turbulent diffusion on the chemistry at the interface between a cold neutral medium (CNM) cloudlet and 
the warm neutral medium (WNM). 

Methods. We perform moving grid, multifluid, ID, hydrodynamical simulations with chemistry including thermal and chemical dif- 
fusion. The diffusion coefficients are enhanced to account for turbulent diffusion. We post-process the steady-states of our simulations 
with a crude model of radiative transfer to compute line profiles. 

Results. Turbulent diffusion spreads out the transition region between the CNM and the WNM. We find that the CNM slightly expands 
and heats up: its CH and H2 content decreases due to the lower density. The change of physical conditions and diffusive transport 
increase the content in the CNM which results in increased OH and H2O. Diffusion transports some CO out of the CNM. It also 
brings H2 into contact with the warm gas with enhanced production of CH^, Hj , OH and H2O at the interface. O lines are sensitive 
to the spread of the thermal profile in the intermediate region between the CNM and the WNM. Enhanced molecular content at the 
interface of the cloud broadens the molecular line profiles and helps exciting transitions of intermediate energy. The relative molecular 
yields are found higher for bigger clouds. 

Conclusions. Turbulent diffusion can be the source of additional molecular production and should be included in chemical models 
of the interstellar medium (ISM). It also is a good candidate for the interpretation of observational problems such as warm H2, CH* 
formation and presence of Hj . 

Key words. Astrochemistry - Turbulence - Diffusion - ISMxlouds - ISM: molecules - ISM: lines and bands 
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1. Introduction 

In the thermally bistable neutral ISM, small structures of cold 
neutral gas (CNM for cold neutral medium, T a; 70 K and 
n ^50 cm"^ ) are immersed in the warm neutral medium (WNM, 
T 8000K and n » 0.4 cm ■^). Both phases are found to have ap- 
proximately the same pressure (iHeiles & Trolandll2003l) . In this 
paper, we investigate how turbulent diffusion processes can af- 
fect the chemistry and the structure of an interface between a 
CNM cloudlet and the WNM. 

Much of the Icnowledge on the cold diffuse medium has been 
established from the detection of molecules, radicals and molec- 
ular ions. Despite its relatively low density the diffuse interstel- 
lar medium is chemical ly active and furthermo re, its chemistry 
is not fully understood dSnow & McCal]ir2006l and references 
therein). The presence of the reactive ions CH"*^ and H:^ together 
with radicals like CH and C2 and molecules like CO in the dif- 
fuse ISM can not be understood with s teady state chemical mod- 
els. For example. iLe Petit et al.l(l2004l) present a detailed model- 
ing of the ISM along the line of sight of ( Per which requires 
three components along the line of sight : a long (4 pc) diffuse 
region at 60 K, a tiny (100 AU) dense region at 20 K, both with 
a large ionization rate, and a series of shocks which account for 
the abundance of CH* and the excited rotational populations of 
H2. This example is not unique and demonstrates the necessity of 
coupling self-consistent models of the structure of the ISM with 
the chemistry. It also shows that CH* and H3 are key species for 
understanding the coupling of turbulence, dynamics and chem- 
istry in the diffuse ISM. 

The ubiquitous presence of CH* and excited rotational lines 
of H2, which are not expected at thermal equilibrium given 



the average temperatures of the cold neutral medium, indicates 
that non equilibrium processes must also be taken into account, 
which will bring enough energy locally to trigger the formation 
of CH* and could contribute to the excitation of H2. is mainly 
formed with H2 and H^ (itself resulting from the ionisation of 
H2 by a cosmic ray) and efficiently destroye d by dissociative 
recom bination with electrons. As reviewed bv ISnow & McCall 
(I2OO6I) . diffuse cloud models do not predict H^ to be abundant, 
except when the ionisation rate is high enough to maintain a 
large formation rate. Steady state chemistry predicts that the 

density scales as: n(ll^) =y^^- Because the H^ density is ap- 
proximately constant and proportional to f , the presence of large 
H3 column densities requires long pathlengths through the CNM 
(4 pc for ( Per as modelled by Le Petit et al. 2004 and using 
^ = 25 X 10"'^ i"') which significantly exceed the expected size 
of such structures (e.g. around 0.1 pc up to 1 pc in the simula- 
tions bv I Audit & Hennebellell20(j5l) . The problem with the path 
length is even more severe because the determined column den- 
sity of H3 towards ^ Per requires a path length larger than 60 pc 
when using the standard value for the cosmic ray ionisation rate 
, ^ ~ 3 X 10"'^ i"' ! As CNM cloudlets are never isolated but 
always immersed in the WNM, we investigate in this paper the 
effect of turbulent diffusion at their interface as an alternative so- 
lution for enhancing molecular abundances in the ISM. A subse- 
quent paper will present the effects of shocks. 

The chemical importance of turbulent diffusion has been ac- 
knowledged in the context of protostellar cores since the work of 
Bol and & de Jond (11982 ) who realised how extra mixing could 
limit element depletion on grains. Since then, many authors have 
computed the effect of turbul ent mixing on the chemistry of pro- 
tostellar and prestellar cores. [Chieze & Pineau des Foretsl(ll989l) 
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and IChieze et alj (1199 Ih have used a simple two-zone model 
and a mixing prescription between the envelope and the cores 
of these objects. Xie et al. (1991) and |Xie et al. (1995) com- 
puted the isothermal evolution of chemist ry with an enhanced 
diffusion coefficient. IWillacv et al.l (12002!) ' later extended this 
work including H2 photo-dis sociation and gas-grain interactions. 
iRawlings & HartquistI (11997 ^) did similar computations but using 
a uniform pressure and a diffusion term based on numerical den- 
sities rather than a Fick law based on mass fraction. 

In the present study, we transpose this series of work to an- 
other context with physical conditions at lower densities and 
lower visual extinctions, more appropriate to the diffuse neutral 
ISM. We also improve on the previously used numerical meth- 
ods. First, we compute the hydrodynamical evolution of a mul- 
tifluid gas and we treat the gas cooling in a self-consistent man- 
ner, completely coupled to the chemical evolution. Second, we 
model the turbulent diffusion of temperature along with chemi- 
cal mixing: this will prove useful to extend the region where the 
formation of molecules is favoured. 

In section 2 we describe our numerical method. We present 
our results in section 3 and discuss them in section 4. Section 5 
draws our conclusions. 



The total number of reactions is now 138. The full list o f re- 
actions can be found in the last appendix of Lesaffre' (^2002") and 
takes its main roots in the early work of Chieze et al. ( 1998): it 
is designed to account for the abundances of the main cooling 
agents of the ISM such as C+, C, O, H2, OH, CO and H2O. The 
dependence of each photo-reaction rate on the local radiation 
field is modeled with a parameter fi. The local rate is obtained 
by multiplying the rate for the ISRF by the factor exp(-y6A,,). In 
most cases, this simple rule is sufficient to get an accurate model 
for the extinction of the radiation field in the wavelength range 
relevant to a given photo-reaction. 

However, the photo-dissociation of the molecule H2 requires 
a special treatment because of its sel f-shielding properties. For 
this purpose, we use equation (37) of iDraine & Bertoldil (Il996h 
with a constant Doppler parameter of =0.1 corresponding to 
a 0.1 km.s"' Doppler broadening (this small value maximises the 
effect of self-shielding). Note that due to the small total column- 
density across our computational box, we do not need a similar 
treatment for the photo-dissociation of CO. 

In the present simulations we use 1.5 10"^ cm^-'s"' at 300 K 
for the dissociative recombination rate of H^ with a cosmic ray 
ionisation rate of ^ = 5 10"'^ s"'. 



2. Method 

We use the multifluid monodimensional (plane paralle l) mag- 
neto hydrodynamical (MHD) code presented in Lesaff reet alj 
(1200 4). This code fully couples the hydrodynamics and the 
chemistry. It also m akes use of a moving grid algorithm 
dPorfi & Drurvl Il987b which allows to spatially resolve the 
shock discontinuities and the diffusion fronts. In order to model 
the diffuse molecular gas we had to supplement the code with 
the following new processes. 



2.3. Diffusion 

The structure and evolution of a front between a CNM cloudlet 
surrounded by the WNM are mainly determined by diffusion 
processes. Diffusion is due to random motions which transport 
thermal energy or chemical species. It can be caused either by 
microscopic collisions between molecules or by random macro- 
scopic motions of turbulent eddies. 

In the following, we evaluate diffusion coefficients for these 
two modes of diffusion. 



2.1. Radiative transfer 

A radiation field equal to the standard interstellar radiation field 
(ISRF) illuminates the right hand (outer) side of the computa- 
tional box. The extinction at the outer side of the box is very 
small: Aio - 0.01, relevant to WNM conditions. The extinction 
at any point inside the computational box is given by 

A,(r)^ A,o + 1. 5x10-^^ cm^ j Nndl (1) 

where r is the position inside the computational box, R is the po- 
sition at the outer (right hand) edge of the box and Nh is the nu- 
merical density of hydrogen nuclei (hence A^, actually increases 
from right to left). In most simulations of this paper, the total 
Ay(0) across the box is 0.14 mag for a total length of the box 
R = 5 pc: we refer to them as "Small cloud" (or standard) sim- 
ulations. In some other simulations, we used A,,(0) =1.4 mag 
and R = 10 pc ("Big cloud" simulations). The local value of A^, 
is used to determine the heating rate due to the photo-electric 
effect on grains and various photo-chemical rates (see below). 

The dense cloud will therefore naturally be formed at the left 
hand side of the computational box, where the extinction is the 
highest. The WNM will be found at the right hand side. 

2.2. Photo-chemistry 

The reaction network in lLesaffre et al.l ( |2004|) has been extended 
with photo-chemical reactions of dissociation and ionisation. 



2.3.1. Molecular diffusion 

It is a difficult task, still heavily discussed, to determine the dif- 
fusion coefficients in a multi-component plasma. Here, we adopt 
a simple view consistent with our modeling of three different flu- 
ids: neutrals, ions and electrons. 

A typical molecular flight takes place at the thermal velocity 
over a length scale equal to the mean free path. We hence adopt 
diffusion coefficients of the form 

dj — CjA (2) 



with j =n,i,e spanning the neutral, ion and electron fluids, where 
A - l/crN is the mean free path (cr = 10"''' cm^ and is the 
total number density of particles per cm^) and 




is the sound speed in fluid j (k is the Boltzmann constant, Tj is 
the temperature of fluid j and fij its average molecular weight). 
We assume that electrons are dynamically coupled to the ions 
and hence take de := c/,. According to the relative average molec- 
ular weights, diffusion is fastest in the neutral fluid and slowest 
in the ion and electron fluid. It is a bit artificial to adopt the 
same mean-free path for all three components independent of 
energy, but as long as these length scales are much smaller than 
the turbulent diffusion lengths scales, our results should not be 
changed. 



p. LesafFre et al.: Effects of turbulent diffusion on the chemistry of diffuse clouds 



3 



d_ 
dr 



The thermal diffusion term takes the form 
~ dkT, 
dr 



(4) 



in the evolution equation for the energy density for each fluid 
j. Nj is the number density of particles in fluid j. Njdj is Njdj 
interpolated linearly at each interface of the computational zones 
and diffusion fluxes are assumed to vanish at the computational 
boundaries. 

We adopt a similar treatment for the chemical diffusion term 



(5) 



in the evolution equation for the number density A^, of species / 
in fluid j, with Xj the mass fraction of species / and p the total 
mass density. Note that we do not take into account the indi- 
vidual thermal velocities for each species but we rather take the 
average sound speed in each corresponding fluid. This may be 
improved in future work. 

2.3.2. Turbulent diffusion 

Consider the evolution equation of the number density of a 
given chemical species (of mass fraction X) with velocity U, net 
chemical production rate R and diffusion coefficient d: 

dN 

+ V.(NV)^R + V.[dpVX]. 



dt 



(6) 



Our purpose is to find an evolution equation for the quantity 
averaged over a given length scale I, which we denote A^. The 
average of equation (|6]l yields 



ot 



R + V .[dpVX] - V .{nu) 



(7) 



where u = U - U and n - N - N. 

The challenge of turbulent models resides in the difficulty of 
finding a simple but accurate model for the behaviour of high or- 
der moments such as V.(nu). We now sketch a brief justification 
for the crude model we adopted. We compute the time evolu- 
tion of n under the first order smoothing approximation (FOSA) 
which allows us to drop all orders higher than two. The differ- 
ence between equations (|6]l and ^ yields 



|^+V.(„U): 

Ot 



-V.{Nu) + (R-R) 



(8) 



where we neglected second order terms and safely dropped the 
term proportional to d since it is generally small compared to the 
turbulent diffusion. The term R - R arises because of di fferential 
reactivity in turbulent eddies (see iLesaffre et alj|2005h and can 
be neglected provided the turbulent mixing time scales are much 
shorter than the chemical time scales. That is hard to justify, and 
can in fact be wrong for most of the rapidly reacting species, 
but we will neglect this term in the present work for the sake of 
simplicity. Note that we also make the convenient approximation 
R{N) = ^(A^) which is consistent with our neglect of differential 
reactivity. 

We now assume that n is given by the integration of the 
source term in equation (|8]l over a given correlation time r, suf- 
ficiently short to neglect the time evolution of spatially averaged 
values. This leads us to the expression for the turbulent term in 
equation 



Mass conservation of the averaged quantities requires A^ = 
p X. We use this relation to expand the previous relation: 



- V.(nu) = V.[p (tu ® u).VX] + V.[X (tu ® V).(p u)] 



(10) 



V.(«u) = V.[tV.(A^u)u]. 



(9) 



We retain only the first term which is analogous to a diffusion 
term with a tensor turbulent diffu sion coefficient K - ru ® u 
which simplifies to a scalar K - ^tu^ for isotropic turbulence. 
The second term is an advection term that is present only when 
turbulence carries a non zero net mass flux. 

Under these simplifying assumptions, the effect of an under- 
lying turbulent velocity field can hence be reduced to an effective 
diffusion coefficient K. A good disc ussion of t he obs ervational 
uncertainties on K can be found in IXie et all (Il995h . In brief, 
we have a good idea of the rms fluctuations of the velocity field 
from the broadening of spectral lines, but we only have a poor 
understanding of the correlation length scale (or equivalently of 
the correlation time scale) of the turbulence. 

Given these uncertainties, we investigated several possibili- 
ties for the value of K. For all computations with turbulent dif- 
fusion presented in this paper, we keep the same formulation as 
for the molecular diffusion coefficients dj above, but we use a 
fixed correlation length L instead of A and the sound speed c„ in 
the neutrals for all fluids. This implicitly assumes that the turbu- 
lent diffusion coefficient scales as the local average sound speed, 
which may be questionable. On the other hand, it yields smaller 
diffusion coefficients in colder parts of the cloud and it implies 
that the ratio of the thermal pressure to the turbulent pressure is 
a fixed constant, which seems reasonable (in particular, it gives a 
uniform turbulent pressure field when the thermal pressure itself 
is uniform). We also tested the more s imple uniform p rescrip- 
tion K = LxQ.2 km.s ' (similar to IXie et all [19951) but we 
did not find any significant qualitative difference with the results 
presented here. In all our tests, we adopted the same diffusion 
coefficient for the chemistry and the temperature. 



3. Results 

3.1. Numerical setup 

In this section, we investigate the steady state of condensation 
fronts by means of time-dependent simulations evolved over 
very long times. We prepare a cold state (CNM) and a hot state 
(WNM) of thermochemical equilibrium at a given total pres- 
sure Po = 3 X lO-' K.cm^^* for a given extinction Am - 0.01. 
More precisely, we evolve a two-zones isobaric model until both 
zones achieve a steady state of thermal and chemical equilib- 
rium at A,,o - 0.01, hence the initial molecular content is small. 
Our adopted values for the abundance of C and O gas phase are 
1.3 10""* and 4.3 10""* respectively. We then select a length I 
{I = 0.5 pc in for our small cloud simulations and Z = 5 pc 
for our big cloud simulations) for the cold "cloud" which will be 
located to the left (inner) side of the computational box and the 
rest of the box (of total length R=5 pc or 10 pc in the small or 
big cloud simulations) is filled with the hot state. We then let the 
gas evolve chemically and hydrodynamically in the box. This 
procedure essentially sets only the total mass in the box, as the 
size of the cloud and the pressure will later relax according to 
the hydrodynamical, thermal and chemical evolution. We later 
refer to the small cloud simulation without turbulent diffusion as 
our reference run. 
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3.2. Evolution of the small cloud without turbulent diffusion 
(reference run) 

The chemistry first reacts to the new local Ay > A^o, with 
molecules reforming in the cloud. Another transient effect is due 
to the partial pressure of ions and electrons which differ in the 
cold and hot phase due to fractional abundance differences. The 
gradient of partial pressures generates a flow of ions across the 
surface of the cloud. This process is responsible for a transient 
chemistry which we do not describe here, as it mainly results 
from our arbitrary initial conditions. Microscopic (molecular) 
thermal and chemical diffusion then build up a smooth spatial 
transition between the cold and hot phases. 

We give the time evolution of the total column densities of 
several species in figure[T]as an illustration of the very long times 
required to reach a steady state. One after the other, chemical 
reactions reach their steady state according to their own time 
scales. After a time of about 0.1 Gyr chemical equilibrium is 
achieved and the shape of the diffusion front (figure [3]) reaches 
an approximate steady-state. However, it is not stationary yet 
as it is moving from the hot towards the cold medium, with the 
cold medium reaching colder temperatures and higher molecular 
fractions. The cloud shrinks in size as the condensation front 
moves and its density gets higher while the hot medium gets 
hotter and more diffuse. Note that the overall pressure rises while 
the front tries to adjust the total cooling rate to the total heating 
rate, on its way to find a steady state. 

This simulation did not show any sign of reaching a steady 
state after more than 1 Gyr. However, extrapolating the 0. 1 Gyr 
time scale to reach steady state in our lowest turbulent diffusion 
coefficient to 20 times lower diffusion would give a time scale 
of 2 Gyr. This computation may hence not have been carried out 
over a long enough time to probe if this front will eventually 
reach a steady state or not. 

The multifluid nature of our computations is not essential to 
the present study. Nevertheless, it allows us to uncover a poten- 
tially interesting partial pressure effect. A slight neutral/charges 
velocity drift of the order of 0.01 km.s ' persists in the hot phase 
throughout the simulation. This drift is sustained by the partial 
pressure equilibration. Indeed, photo-ionisation is more effective 
in the hot phase and a partial pressure gradient drives a flow of 
charged ions towards the cold medium. This is the cause of a 
current of charges from the hot toward the cold medium : 
ions accumulate at the surface of the cloud where their radiative 
cooling produces a dip in the temperature profile (see figure [3]). 
If self-gravity of the cloud was taken into account, this would 
produce a convective instability at the surface of the cloud that 
would soon smooth out the temperature and profiles. This 

effect might also be responsible for the failure of the front to 
reach a completely static state. 

3.3. Effect of turbulent diffusion 

Although we display the results for big cloud simulations, we 
restrict our analysis to small cloud simulations in this section. 
The big cloud simulation was run following the same numerical 
experiment as for the small cloud, but we investigated only two 
different values for the turbulent diffusion coefficient. Big cloud 
simulations are discussed in section |4n 

We restart our computation from f = 0. 1 Gyr in the reference 
run and switch on turbulent diffusion while resetting the time 
to 0. After about another 0. 1 Gyr has occurred, a static state is 
reached. The time to reach a steady state is much shorter than in 
the case without diffusion, because the effective diffusion time 



scales are much smaller. The ion current is still present, but the 
diffusion smoothes out the dip in temperature and in compo- 
sition. This seems to allow the front to find a final equilibrium. 
The pressure in this final equilibrium turns out to be nearly 40% 
less than in the case without turbulent diffusion at time 0.1 Gyr. 

We now probe the dependence of this steady state on the 
magnitude of the turbulent diffusion coefficient. We compare the 
state of the simulation box at time 0. 1 Gyr for different turbulent 
length scales L ranging from 1.25 x 10-*R to 1.6 x IQ-^R with 
a factor 2 increment. With R = 5 pc, this range extends from 
L - 6 10"* pc to L = 0.08 pc. Note that the mean free path is at 
a scale of around 3 x 10"^ pc in the cold medium, which coiTe- 
sponds to 6 X lO^^R and makes molecular diffusion negligible 
by at least a factor 20 with respect to turbulent diffusion in all 
our computations. The upper bound L = 0.08 pc is roughly a 
seventh of the size of the cold cloud at this level of diffusion (re- 
member L is set as a fraction of the size of the computational box 
but it is difficult to predict what the actual size of the cloud will 
be befor e actually runnin g the simulation). As a point of com- 
parison, IXie et al.1 (1 19951) built their diffusion coefficients with 
L = 0.1 - 0.5 pc and a fixed velocity of 1 km.s ' for a molecular 
cloud of size 1 pc. 




Time (yr) 



Fig. 1. Reference run (small cloud without turbulent diffusion): 
time evolution of H2 (solid, value scaled by a factor 10"*), 
CO (soHd), OH (dotted), H2O (dashed), H+ (dash-dotted), CH+ 
(dash-dot-dot-dot) and CH (long-dash) column-densities (cm"^) 
integrated over the computational box. 



Figure|2]displays the total column-densities of the molecular 
species of main interest for various diffusion coefficients. Except 
for the column-densities of CH which decreases, and those of 
CO and Hs^ which stay at the same level, the effect of turbu- 
lent diffusion is to boost the production of molecular species. 
In the following, we proceed to a more detailed investigation 
of the structure of the front: for each molecule of interest in 
each simulation, we computed which source and sink terms were 
dominant, including reaction and diffusion effects. This reveals 
that H2 molecules mixed into hot gas are the main cause for the 
diffusion-enhanced production of molecules. 
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Fig. 2. CO (solid), OH (dotted), H2O (dashed), H+ (dash-dotted), 
CH^ (dash-dot-dot-dot) and CH (long-dash) column-densities 
(cm~^) integrated over the computational box for condensation 
fronts at time 0.1 Gyr with various diffusion lengths L (the 
big cloud run with the highest L is actually displayed at time 
0.02 Gyr, when the WNM has disappeared). Plus signs are for 
the small cloud and diamonds stand for the big cloud runs. Note 
that the values for the runs without turbulent diffusion at time 
0.1 Gyr have been arbitrarily placed at L = 3 x 10"^ pc (small 
cloud) and L = 0.1 pc (big cloud). 



3.3.1 . Thermal and density structure 

As illustrated in figure [3] the main effect of increasing the diffu- 
sion is to spread out the length of the front. As a result of the heat 
exchange between the hot and the cold phase, the cold phase be- 
comes slightly hotter and the hot phase becomes slightly colder 

As already mentioned, the total pressure at steady state de- 
creases when diffusion is enhanced. Besides, the cold medium 
temperature increases slighlty in the case with diffusion. As a 
result, the total density in the cloud decreases by more than a 
factor 2 in the case with diffusion (see figure |4|. The cloud ex- 
pands to keep its total column-density constant (remember that 
the total mass of the simulation box remains constant). 

Everywhere in all simulations the dominant heating process 
is the photo-electric effect. The main cooling processes are from 
the inside of the cloud to the WNM (left to right): cooling, 
diffusive cooling (on the WNM side of the diffusion front) and 
H (Lyman a) cooling. In the WNM, O cooling is also impor- 
tant, but amounts to at most a half of the Lyman a cooling. 
When turbulent diffusion is present, molecules are abundant in 
the WNM (see below) and H2 cooling sets in in between the re- 
gion where diffusive cooling and Lyman a cooling prevail. Note 
that although O collisional excitation b y electr ons is adopted, we 
use H2 cooling functions from lLe Bourlot et al.l (1 19991) which in- 
clude excitation by H, He and H2 only. 

3.3.2. H2 abundance 

The abundance of H2 usually results from the balance between 
its formation on grains and its destruction by photodissociation. 
However, for cases with turbulent diffusion, in the intermediate 
region between the CNM and the WNM, the main source of H2 
is the divergence of its diffusive flux. 




0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 
Distance (pc) 



(b) Tn (Big cloud) 
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Fig. 3. Temperature (K) structure of the front at steady-state for 
various diffusion coefficients: (a) upper panel (standard run) the 
solid line is for molecular diffusion, the dotted line for L - 6 x 
10 pc, the dashed line for L = 10"^ pc and the dash-dotted line 
for L = 8 X 10"^ pc ; (b) lower panel (big cloud run) the solid 
line is for molecular diffusion, the dotted line for L - 0.16 pc 
and the dashed line for L = 0.8 pc. 



An increase in the diffusion coefficient considerably en- 
hances the abundance of H2 in the diffusion front (see figure 
lU. This is mainly a diffusion effect: diffusion mixes the cold 
molecular gas from the cloud with the external hot atomic gas. 
The diffusion flux provides a source of H2 molecules in the hot 
phase which is strong enough to counteract photo-dissociation 
and manages to sustain high concentration levels of H2 in a hot 
medium. Hence turbulence brings H2 molecules in contact with 
hotter gas: this effect will be the source of a lot of efficient molec- 
ular chemistry. 

The decrease in H2 abundance in the cloud for enhanced dif- 
fusion results from the combined effects of the general spread 
in H2 (diffusion pumps H2 from the cold cloud) and the lower 
density (H2 formation on grains is less efficient). 
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(b) nH (Big cloud) 
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Fig. 4. Density (cm ^) structure of the front for various diffusion Fig. 5. Abundance profile of H2 for various diffusion coefficients 



coefficients (same labels as figure [3]l. 



(same labels as figure |3]l. 



3.3.3. CH+ abundance 

The most spectacular consequence of increasing the tempera- 
ture of H2 rich gas is the production of the CH^ molecular 
ion. In absence of diffusion, CH+ is produced via the reaction 

+ H — > CH^ and destroyed by photodissociation in the hot 
gas or by CH^ + H — > + H2 in the cold gas. 

However, the reaction + H2 ^ CH^ + H can be a much 
more efficient source for CH^ provided that the 4640 K activa- 
tion temperature is overcome in presence of H2. As clearly seen 
in figure|6] diffusion provides an excellent means of bringing to- 
gether high temperature and H2 rich gas in order to make CH^ 
through this channel. 



3.3.4. H3+ abundance 

Hs^ is made via the reaction H2^ -H H2 — > Hj,^ + H and de- 
stroyed through dissociative recombination with electrons. As a 
result, n{H^) = (n{H2)/ki,n(e) remains valid throughout most 
of the region where is present. The recombination rate kg is 

a decreasing power of the temperature (it behaves like 1/ y/T). 
Furthermore, it turns out that the electron abundance n(e) also 
scales approximately like in the intermediate region. One 



should therefore expect more in the hot phase. However, 
photo-dissociation usually destroys the source of H2 molecules 
in the hot phase. Turbulent diffusion allows to replenish part of 
the high temperature gas in H2 molecules. This effect produces 
a peak of Hs^ inside the diffusion front as illustrated in figure |7] 

3.3.5. OH abundance 

In our simulations, OH usually results from the competition be- 
tween the recombination of H20^ and photo-dissociation (for 
the conditions where the H2 fraction is greater than about 0.2- 
0.3, H20^ is destroyed by H2 and recombination of H^O^ be- 
comes the main source, Pineau des Forets 2007, private commu- 
nication). However, the H-atom transfer O -1- H2 — > OH + H can 
become the dominant source if its 2980 K activation temperature 
is overcome in presence of H2. Figure[8]shows that this is indeed 
the case at the CNM/WNM interface when turbulent diffusion is 
able to bring enough H2 at sufficiently high temperature. 

An other noticeable effect of turbulent diffusion is a fac- 
tor 3 to 5 rise of the relative abundance of OH in the cold 
cloud from the case without turbulent diffusion to the case with 
turbulent diffusion. This is less straightforward to understand. 
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Fig. 6. Abundance profile of CH'*' for various diffusion coefii- Fig. 7. Abundance profile of for various diffusion coeffi- 
cients (same labels as figure [3]l. cients (same labels as figure |3]l. 



Indeed diffusion slightly increases the temperature and strongly 
decreases the density in the cold cloud which is expected to 
decrease the recombination rate of H20^, hence lowering the 
OH abundance. However, H20^ comes from H^ via the reac- 
tion H+ -H O ^ 0+ + H followed by OH+ + H2 ^ H2O+ + H. 
The former reaction is subject to a 227 K activation temperature, 
which rate increases by at least a factor 6 in the cold cloud be- 
tween the case with and without diffusion. H^ is itself subject to 
at least a factor 3 rise in relative abundance in the cold medium 
due to its lower recombination rate thanks to a slightly higher 
temperature and a less dense gas (the electron fraction remains 
roughly the same, as it is bound to the abundance of C^). The 
net effect is a rise in OH relative abundance due to the change in 
the physical conditions inside the CNM. 

The very sharp dip in OH relative abundance at the surface 
of the cloud when turbulent diffusion is switched off (see the 
solid line on figure|8]l is in fact composed of two dips. The inner 
one is caused by the dip in temperature and the sensitivity of the 
H^ H- O reaction rate to temperature. The outer one is caused by 
the decrease of H2 that lowers the H20^ rate of production in 
a region where O H- H2 is not activated yet. Only the outer dip 
remains for L = 6 x 10"^ pc (see the dotted line on figurelS). 



3.3.6. CO abundance 

In the densest parts of the cloud, CO is produced through the 
reaction HCO^ H- e — > CO + H. That reaction accounts for the 
inner 0. 1 pc in the small cloud cases and the inner 3 pc for the 
big cloud with the highest turbulent diffusion coefficient (in the 
latter case, CH -i-O — > CO + H also prevails in the inner 0.5 pc). 
HCO^ is indeed pointed out as one of the most likely sourc es for 
CO in diffuse and translucent clouds (' Liszt & Lucas|[T998h . 

In the rest of the simulation box, the molecule CO is most of- 
ten produced by C^ + OH CO -I- H'*" and destroyed by photo- 
dissociation. In these regions its relative abundance is hence di- 
rectly linked to the relative abundance of OH and density, as 
demonstrated by figure |9] In a small region at the right hand side 
of the diffusion front, the reaction CO^ + H -> CO + H^ is the 
most efficient source of CO. For high diffusion coefficients, the 
abundance profile of CO can also be subject to diffusive mix- 
ing that removes CO from the CNM. The enhanced CO in the 
front, due to enhanced OH, is never strong enough to produce a 
significant local maximum in the diffusion front. 
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(b) OH (Big cloud) 
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Fig. 8. Abundance profile of OH for various diffiision coeffi- Fig. 9. Abundance profile of CO for various difiiision coeffi- 
cients (same labels as figure [3]i. cients (same labels as figure |3]l. 



3.3.7. H2O abundance 

H2O is generally produced by the recombination of the molec- 
ular ion H30^, which comes from proton exchange between 
HaO^ and H2 and is hence linked to the abundance of as 
is OH in the cold cloud. H2O experiences the same rise in abun- 
dance as OH in the cold cloud (see figure [TOli but has an extra de- 
pendence on the H2 fraction. This effect might be responsible for 
the actual decrease of the H2O relative abundance in the CNM 
at high diff'usion coefficients. However, the dominant source of 
H2O in the diffusion front is proton exchange between H2 and 
OH, as both molecules are abundant in the front for the reasons 
stated above. 

3.3.8. CH abundance 

CH results from the dissociative recombination of CH3^ and is 
mainly destroyed by photo-dissociation. CH3^ is the result of 
two successive proton exchanges with H2 starting from CH^. 

The decrease of CH in the cloud is related to the decrease 
of the H2 abundance (see figures |5] and fTTT i. In the diffusion 
front, the enhanced CH^ abundance does not compensate for the 



higher temperature (which lowers the recombination rate) and 
the lower abundance of H2. 

The activation temperature of 1.41 10'* K for the reaction 
C H- H2 — > CH -H H is too high to be overcome in the same way 
as for CH^. Hence there is a net decrease of CH column-density 
when the diffusion is enhanced (the decrease comes from the 
lower H2 abundance in the cloud). 

3.3.9. Line profiles 

We post-processed our simulations using a crude radiative trans- 
fer code in order to probe the effect of diffusion on the line pro- 
files. We assume that a given species in each parcel of gas in the 
simulation emits (or absorbs) light with a Gaussian line shape 
in velocity space. The mean of the Gaussian is set to the local 
speed of the fluid, its width is set to the thermal speed of the 
species and its amplitude is set according to the column density 
of the parcel. Note that we implicitly assumed that the turbulent 
velocity is a constant fraction of the local sound speed. In this 
simple picture such turbulent motions do not modify the shape 
of the spectra if the velocity scale is expanded appropriately. We 
hence neglect the turbulent Dopper shifts and concentrate only 
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(a) H20 (Small cloud) 
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Fig. 10. Abundance profile of H2O for various diffusion coefii- Fig. 11. Abundance profile of CH for various diffusion coeffi- 
cients (same labels as figure [3]i. cients (same labels as figure [3]l. 



on the intrinsic line shapes which reflect thermal and chemical 
changes. 

In a first step, we don't take into account excitation effects. 
Assuming an optically thin medium, or a line profile seen in ab- 
sorption, we compute the line shape across the box as the total 
column-density of species M per unit velocity bin centred on 
velocity v: 



NAM 



» = r dr N{M) exp [- J-* ] 
Jo 



2kr 



27rkr 



(11) 



where N(M) is the number density of species M, fi is its molec- 
ular weight, u is the local velocity of the fluid associated to M 
and T its temperature. Note that the integral of Nv{M, v) from 
V = -00 to V - +CO yields the total column-density of species M 
across the simulation box. Taking turbulence into account could 
be done simply by artificially decreasing 

Figures[T2]and[T3]display A^,, for a few molecules and atoms 
in the case without turbulent diffusion and in a case with turbu- 
lent diffusion. Without turbulent diffusion, all molecules show 
a narrow Gaussian profile corresponding to the temperature of 
the CNM. Indeed, molecules are absent from the WNM. 
and O on the contrary are equally abundant in both phases and 



they show a wider secondary Gaussian feature, corresponding to 
emission from the WNM. The amplitude of this second Gaussian 
is smaller because there is much less total column-density in the 
WNM than in the CNM. 

When turbulent diffusion is present, molecules are formed 
in the transition layer at intermediate temperatures between the 
CNM and the WNM. This enhances the contribution of warm 
gas in the integral and gives rise to broad line wings. The ef- 
fect is strongest in the case of CH^ which has the most effi- 
cient yield due to turbulent diffusion. The effect is inexistent 
for CH which subsists in the cold phase only. This provides a 
potentially powerful comparative probe between CH^ and CH. 
Indeed, C H^ lines are bro ader than CH and CO lines in the dif- 
fuse ISM (iPan et alJl2004l) . 

We now compute the emission profiles of a few atomic tran- 
sitions for which our code provides the level populations at 
steady-state (they are used to compute the atomic cooling func- 
tions). We integrate equation ( fTTT i using the local emissivities of 
these transitions instead of the quantity N{M). The results are 
displayed in figure [141 The fine-structure transitions of O show 
a strong effect with nearly exponential wings appearing in the 
case with turbulent diffusion. O has a nearly constant abundance 
across the simulation box. The wing effect in this case is en- 
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Fig. 12. Line profile for H2O, CH+, CO and CH (see text). The 
solid line is in the case without turbulent diffusion, the dashed 
line is in the case with turbulent diffusion (L = 8 x 10"^ pc). 
A Gaussian line has been adjusted on the core of each profile 
(for velocities lower then 0.2 km.s '). It is displayed as a dotted 
line for the case without diffusion and as a dash-dotted line for 
the case with diffusion. In most cases, the dotted line cannot be 
distinguished from the solid line. 



tirely due to the structure in density and temperature. The fine 
structure transitions have energy temperatures of 227 K, 98 K 
and 326 K intermediate between the CNM and WNM tempera- 
tures. Besides, the critical density for O is much larger than the 
CNM density. It turns out that the emissivity of these O lines is 
highest in the intermediate temperature and density region inside 
the front. The change of integrated emissivity in the line hence 
reflects only the change in total column-density found in this in- 
termediate region only. O line profiles are a good probe of the 
thermal structure of such condensation fronts broadened by tur- 
bulent transport of heat. On the other hand, the 158 fim (92 K) 
transition of C^ has a much lower critical density and its inten- 
sity is highest in the CNM : its overall line profile barely sees the 
structure of the front. 



4. Discussion 

4.1. Initial conditions 

Our assumption of a closed box model determines the total mass 
contained in our simulations. In particular the column-density of 
hydrogen nuclei is determined by our initial conditions (ie: 8 x 
10''^ cm"^ in the standard runs and 8x10^" cm"^ for the big cloud 
runs). Then, the properties of the radiation field at the outer side 
of the box completely determine the final steady-state. Because 
our ID simulations are unable to self-consistently generate the 
turbulent motions, the value of the diffusion coefficient also has 
to be treated as an external parameter. 

As has been detailed above, the total molecular enhance- 
ments due to turbulent diffusion can either take place in the 
whole cloud or be located at the interface. The decrease of H2 
(and hence CH) is caused by the lowered density of the cold 
cloud which happens over the whole cloud. The enhanced pro- 
duction of CO, OH and H2O due to diffusion of H^ ions from the 
WNM to the CNM occurs over the whole cloud as well. Finally 



Fig. 13. Line profile for H2, Hs^, C^ and O (see text). Labels are 
the same as in figure [T2l 
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Fig. 14. Emissivity profiles for the 63 fim (227 K), 145 //m 
(98 K) and 43 /rni (326 K) fine structure transitions of O and the 
158 fim (92 K) fine-structure transition of C^ (see text). Labels 
are the same as in figure [T2l 



H^ and CH^ production needs high temperature and presence of 
H2 molecules which can happen only at the interface. We there- 
fore expect different scalings for these effects with respect to the 
size of the cloud. 

In order to probe the effect of increasing the size of the cloud 
on these global enhancements, we computed big cloud models 
with 10 times as much mass as in our standard run. We simulta- 
neously increased the turbulent diffusion coefficient to account 
for the bigger size of the cloud: bigger clouds should exist in re- 
gions where the correlation length of the turbulence is larger. We 
investigated a factor 2 enhancement and a factor 10 enhancement 
of the diffusion coefficient compared to the largest value we used 
in the small cloud case. In the latter case, the simulation box of 
10 pc was probably not large enough to contain the wole cloud 
extended by turbulent diffusion and the WNM phase almost dis- 
appears. In that last case, most of the volume surrounding the 
cloud would be occupied by the gas at intermediate temperature. 
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Table [T] summarises the results for some species of interest. The 
molecular enhancements remain qualitatively the same for the 
big cloud. However, the magnitude of the yields greater than 1 
increases when the size of the cloud increases: non-linear effects 
in the molecular enhancement make diffusion more efficient for 
the bigger cloud. 

We also report in table [T]the measured column densities to- 
ward ^ Per Although the models were not aimed at fitting the ^ 
Per data, we obtain a relatively good agreement for most species 
for the Ctio model, which is run with a standard CR ionisation 
rate (5 x 10"'^ s"'). Diffusion considerably enhances the size of 
the region where is present, with a pathlength of 10 pc in 
this case. In particular, we produce reasonable quantities of H3 
(1.5 10'^ cm~^ for our best case) without having to increase the 
cosmic ray ionisation rate as observations sometimes seem to re- 
quest (.McCall et alJl2002t [Le Petit et alJ[2004l: ISnow & McCalll 
I2OO6I) . Our models compared to the ^ Per observations also fall 
short of CO and CH^, both of which probably require the pres- 
ence of an MHD shock. Furthermore, in the case of the big 
cloud, CO photodissociation may need a treatment for its self- 
shielding, not included in the present simulations yet. 

4.2. Observable effects on the lines 

The effect on the global (ie: neglecting excitation effects) line 
profiles is not very strong, but the recent advances in observa- 
tional techniques may allow the detection of subtle variations in 
the line shapes. However, for a typical molecule, rovibrational 
states provide a plethora of possible transitions that will be trig- 
gered at a variety of different temperatures. 

We illustrated this effect for O: its line profile shows no dif- 
fusion effect when excitation is neglected (figure [T3T l but there 
is a clear difference for individual lines (figure [T4ll. Such transi- 
tions should be good observational probes. Furthermore, a slight 
increase in the wings of H2 (see figure [T3T l could mean a huge 
impact for the lines found at intermediate temperatures in the 
front. In fact, turbulent diffusion of H2 towards warm gas could 
be an explanation for the broadening of H2 lines as seen by 
lLacouretal1(l2005h . 

4.3. Warm H2 

Several authors hav e reported observations of warm H2 in cold 
molecular gas (see iGry et al] 120021: iFalgarone et al.ll2005 l for 
example). Turbulent diffusion could potentially be very efficient 
at increasing the column-density of H2 seen at intermediate tem- 
peratures, which would hence be seen as warm H2. An addi- 
tional hint that turbulent diffusion might be relevant in this prob- 
lem is the presence of CH^ in the line of sigh ts with warm H2 
(lFrisch& Juralll980t iLambert & Dankslll986h . Models of pho- 
ton dominated regions (PDR) with turbulent diffusion should 
hence account much bet ter for the observed H2 excitation di- 
agrani s as published by iGry et alJ (l2002h and IFalgarone et al.l 
(I2005h . 

However, the computation of the whole spectrum needs to 
be carried out with a multigroup radiation transfer code in order 
to assess this effect quantitatively. This has been done already in 
the framework of steady-state PDRs computations which on the 
other hand lack a self consistent model for their hydrodynami- 
cal background. Including diffusion in these simulations could 
be one way of proceeding further. Note that diffusion would par- 
tially remove the uncertainty on the level population of H2 newly 



formed on grains since in the case of turbulent diffusion, H2 pro- 
duction is quite often dominated by transport processes. 

4.4. CH^ and CH formation 

Spectroscopic measurements of CH^ and CH usually quote 
column-densities of the order of 10'^ to 10^^ cm"^ for 
each velocity component, with CH and CH^ approxi mately 
correlated Ca mbert et al. 1990; Penprase 1993; Allen 11994 : 



iCrennv & FedermanI l2004t iPan et al.1 12004t I2005D . Such high 
amounts of CH^ were a puzzle for chemical modelers because 
the conditions for efficient CH^ formation were hard to be met 
in molecular gas. Since then, several propositions were made to 
overcome this problem. Drift velocities are an elegant way of in- 
creasing the reaction temperature between charged and neutral 
species without ne cessarily increasing the therm al energy of the 
gas: MHD shocks dPineau des Forets et aT]| 19861) a nd MHD tur- 
bulence dissipation in vortices ( Joulain et al.l 19981) make use of 
this fact. 

The idea that CH^ would be produced in turbulent boundary 
layers at the interface between cold cl ouds and the i r surro unding 
medium was originally suggested bv lDulev et alj (|1992|) . based 
on steady-state chemistry computations. Our simulations are the 
first ones to test this idea. In our model we are unable to account 
for column-densities as high as observed (the most we produce is 
6.2 10" cm - for our big cloud run). However, lines of sight with 
a small amount of CH^ (< 10'^ c m~^) and a normal amount of 
CH (2 10'^ cm-2) are also found (GredeI' 2004l;lPan et al.ll2004h 
which may be best explained by our turbulent mixing scenario. 

CH^ lines are broader than CH fines (ILambert et alj 119901 ; 
ICrawford et aP IT991 iMceetan I200TI; |Pan et all l2005l) . This 
might be a way to discriminate between the various CH^ forma- 
tion channels. Indeed, MHD shocks would produce little ther- 
mal broadening as the drift velocity is the main factor used to 
overcome the activation temperature. Turbulence dissipation will 
produce some line broadening as both heating and drift motions 
in the vortices are needed to produce CH^. Turbulent mixing 
finally will produce the largest broadening as the gas has to be 
readily brought up to temperature higher than the activation tem- 
perature in order to produce CH^ . 

Fi nally, density probes such as the amount and excitation 
of C2 ( lGredellll999i |2004 produced may also be used to dis- 
criminate between the various scenarii as the densities will tend 
to be high in MHD shocks and low in CNM/WNM interfaces. 
However, we need to improve the present work with an extensive 
chemical network in order to check whether C2 is indeed formed 
at the interface in association with CH^ or if it would reside in 
the cloud in our scenario. Note that our model predicts that 
should correlate with CH+ to some degree if both are formed at 
the interfaces of CNM clouds. Unfortunately, the present surveys 
of H3 are too scarce to test this idea. 



4.5. Self-sliielding 

As mentioned earlier, the Doppler parameter that controls the 
e ffect of the broadening of the lines on the self-shielding of H2 
(IDraine & BertoldilfT996 ) was implemented with a small value 
of 0.1 km.s '. We ran another simulation without turbulence 
with a more reasonable value of 1 km.s"' but the resulting to- 
tal steady-state column-densities were at most 7% lower than in 
our reference run. This shows that such a variation of the self- 
shielding function of H2 has a small impact on our results, as 
expected from the high level of saturation of H2 lines. 
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Species 


c 


ft 


>• 


C 


Ct2 


y2 


Ctio 


yio 


^ Per 


H 


4.6(19) 


7.5(19) 


1.62 


1.1(20) 


1.8(20) 


1.66 


1.8(20) 


1.64 


5.7(20) 


H2 


2.0(19) 


6.9(18) 


0.34 


3.6(20) 


3.4(20) 


0.93 


3.5(20) 


0.96 


3.7(20) 


C^^ 


1.6(16) 


1.6(16) 


1.02 


1.4(17) 


1.5(17) 


1.05 


1.5(17) 


1.07 


1.8(17) 


C 


4.3(13) 


1.5(13) 


0.34 


8.1(15) 


4.1(15) 


0.51 


4.6(15) 


0.57 


2.9-3.6(15) 





3.7(16) 


3.8(16) 


1.02 


3.6(17) 


3.6(17) 


1.02 


3.7(17) 


1.05 




CO 


2.7(11) 


2.1(11) 


0.78 


2.3(13) 


9.5(12) 


0.41 


1.5(13) 


0.65 


5.4(14) 


OH 


4.8(11) 


2.1(12) 


4.40 


2.0(12) 


1.1(13) 


5.65 


1.7(13) 


8.53 


4.0(13) 


H2O/OH 


7.9(-2) 


3.6(-2) 




8.0(-2) 


9.8(-2) 




1.3(-1) 








7.5(10) 


9.5(10) 


1.27 


8.7(11) 


3.0(12) 


3.46 


1.5(13) 


17.25 


8.0(13) 


CH+ 


7.2(09) 


9.0(10) 


12.56 


3.6(10) 


1.8(11) 


5.11 


6.2(11) 


17.37 


3.5(12) 


CH 


6.1(11) 


4.1(10) 


0.07 


3.9(13) 


1.6(13) 


0.41 


1.8(13) 


0.46 


1.9-2.0(13) 


H2O 


3.8(10) 


7.6(10) 


2.00 


1.6(11) 


1.1(12) 


6.89 


2.1(12) 


13.33 





Table 1. Total column-densities (in cm"^) of species of interest with their relative yields, c and C are the column-densities across 
the small and big cloud calculations without turbulent diffusion, c, is for the small cloud when maximum turbulent diffusion is 
included with L = 8 x 10 ~ pc. Cti and Ctio are the results for the big cloud and turbulent diffusion with a correlation length L2 - 2L 
and Lio - lOL respectively. The relative yields are computed as y = Ct/c, y2 = Cti/C and yio - Cao/C. Note that the yields are 
not strictly 1 for (H-1-2H2), O and : this gives an estimate of the mass conservation error in the code for individual species. 
Column-densities as observed in the ^ Per line of sight are provided in the last column. The last row diplays t he ratio of H2O to OH 
column-densities. The last column contains the observed column-densities as quoted bv lLe Petit et al.l (|2004|) . 



4.6. Diffusion coefficients 



We already discussed (see section 12.3.1b the uncertainty on the 
turbulent diffusion coefficients. We concentrate here on the pos- 
sible difference between diffusion coefficients for different fluids 
and for heat diffusion. As already mentioned, diffusion coeffi- 
cients in a multicomponent plasma are poorly known. Electrons 
should diffuse more rapidly than the ions due to their much lower 
molecular weight, but they are strongly coupled to the ions due 
to electrostatic forces. As a result, it is likely that the effective 
diffusion coefficient for the ions+electrons fluid is higher than 
the value c/, we adopted. Also, the long range of electrostatic 
forces helps transferring energy without necessarily moving par- 
ticles. As a consequence heat diffusion coefficients are generally 
higher than particle diffusion coefficients. 

On the other hand, turbulent diffusion coefficients depend 
mainly on the properties of the turbulent flow (as hinted at by 
section 12.3. It . As a result, turbulent heat diffusion and parti- 
cle diffusion are likely to proceed with the same coefficients. 
However, MHD turbulence could provide different properties 
for the charged and neutral fluids. We ran simulations with en- 
hanced diffusion in the charged fluid that showed higher yields 
than in the case c/, = de- Indeed, an enhanced diffusion coeffi- 
cient in the charged species introduces additional molecular pro- 
duction in the CNM via H+ transport from the hot to the cold 
medium because OH production is highly se nsitive to the H^ 
ion density. Magnetis ed shocks (Pineau des Forets et al.l [19861) 
and MHD turbulence jJoulain et aL 19981) can also induce drift 
velocities between the charges and the neutrals that may help the 
production of CH^-like species which involve reactions between 
the two fluids and an activation temperature. 



4. 7. Numericai diffusion 

In any non-Lagrangian code, advection errors give rise to some 
numerical diffusion. The coefficient for numerical diffusion can 
be estimated as a fraction of D = VrAr where Vr is the relative 
velocity of the grid with respect to the gas and Ar is the res- 
olution of the grid. In the present study, the resolution is very 
high and grid velocities are very slow because the gas evolves 
over chemical and diffusion time scales: the numerical diffusion 



is always negligible except in the very early stages when the grid 
itself finds its way to the diffusion front. 

Nevertheless we wish to warn hydrodynamical modelers 
willing to include molecular reactions in simulations: they 
should be careful with their numerical diffusion. Indeed, for 
strongly dynamical flows, advection fluxes will be much greater 
than the diffusion fluxes and spurious diffusion much greater 
than the molecular diffusion might occur. In some cases, spu- 
rious molecular production might result. 

5. Conclusions and future prospects 

Turbulent diffusion affects the abundance of OH, H2O, CH^ 
and H3 the strongest. In our most favourable case, the column- 
density increases by more than a factor of 10 compared to the 
steady state model. This shows that some dynamical effects may 
lead to significant abundance variations, of comparable magni- 
tude as uncertainties in reaction rates, and should therefore not 
be ignored in ISM modelling. 

We have shown that H2 mixed into hot gas can be the source 
of warm molecular production thanks to reactions that otherwise 
do not reach their activation temperature: we suggest this sce- 
nario as a new possible channel for the formation of molecules 
such as CH^ or H3 but also OH and H2O in the big cloud case. 
Conversely, the change of physical conditions within the CNM 
favour a higher H^ abundance which can be the seed for OH 
and H2O production (and to a lesser extent, CO). We verified 
that these effects are robust by varying the parameters in our 
simulations such as the total column-density of the cloud or the 
self-shielding function. 

Furthermore, line profiles of transitions with energy differ- 
ences intermediate between CNM and WNM temperatures can 
be considerably affected even for species with no total column- 
density yield (such as O). In particular, models for lines of sight 
containing warm H2 should benefit from a turbulent diffusion 
treatment. In the future, a combination of observational diag- 
nostics should provide good independent probes for the thermal, 
diffusion and chemical processes at play in the transition regions 
between CNM and WNM. For example, some species like O or 
CH are insensitive to interface chemical effects. O excitation is 
sensitive to thermal effects only. Comparative spectroscopy be- 
tween lines of O, CH, CH^, H2O could hence be a very efficient 
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tool to disentangle these different phenomena. A typical signa- 
ture of such diffusion fronts would be the presence of CH^ with 
little CH and wings in the line profiles of CH+, O and H2O but 
not in CH nor C+. 

Thanks to our multifluid treatment, we have also discovered 
a marginal effect due to partial pressure equilibration. The higher 
ionisation degree in the WNM is responsible for a large scale 

ion currents towards the cloud. These could in principle be 
the source for some magnetic field and may also induce thermal 
inversions in the temperature profiles of self-gravitating clouds. 

Finally, we want to insist on the necessity for simulations 
including thermal and chemical effects to account for these tur- 
bulent mixing effects. At a larger s cale, 3D numerical simula- 
tions have already been undertaken dEsquivel et al.ll2006i for ex- 
ample) to probe turbulent mixing effects on the interface be- 
tween the warm ionised medium (WIM) and the hot ionised 
medium (HIM). We wish to motivate multidimensional compu- 
tations of coupled hydrodynamics and chemistry that could help 
cahbrate the turbulent diffusion coefficients and assess the con- 
ditions when a turbulent diffusion formalism is appropriate. In 
the near future, such studies could design subgrid models for 
turbulent mixing in astrophysical objects other than stars. 

In this work, we ran our simulations until they reach a steady 
state, which may take up to 0. 1 Gyr when turbulent diffusion is 
included. This provides a technical point of reference which al- 
lows us to filter all transient chemistry effects and leave only 
the turbulent diffusion effects. However, the ISM is shaken by 
violent events such as supernovae shocks on much shorter time 
scales of the order of 0.01 Gyr. In a forthcoming paper, we plan 
to investigate the effect of hydrodynamical and MHD shocks 
running through the structures studied in this paper. We will 
hence be abl e to put a more consi stent dynamical picture behind 
the model of Le Petit et alj(l2004l) . With the additional chemistry 
coming from compression and neutral/charge drift, we should 
also improve our comparison to the observations. 
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